function [b]=banbks(a,m1,m2,al,indx,b)

% [x]=banbks(a,m1,m2,al,indx,b)
%
% band solver forward and back substitution of a al x = b
%
% a - the upper triangular portion of the matrix
% m1 - number of columns before the diagonal
% m2 - number of columns after the diagonal
% al - the lower triangular portion of the matrix
% indx - the index of rows for reordering
% b - the rhs of the equation
%
% x - the solution to (al a)^-1 b = x

n = size(a,1);

mm=m1+m2+1;
l=m1;
for k=1:n
   %Forward substitution, unscrambling the permuted rows as we go. 
   i=indx(k);
   if(i~=k)
      dum=b(k);
      b(k)=b(i);
      b(i)=dum;
   end
   if(l<n)
      l=l+1;
   end
   for i=k+1:l
      b(i)=b(i)-al(k,i-k)*b(k) ;
   end
end

l=1;
for i=n:-1:1 
   %Backsubstitution.
   dum=b(i);
   for k=2:l
      dum=dum-a(i,k)*b(k+i-1);
   end
   b(i)=dum/a(i,1);
   if(l<mm) 
      l=l+1;
   end
end
